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Motivated by impurity-induced magnetic ordering phenomena in spin-gap materials like TICuCIs, we 
develop a mean-field theory for strongly disordered antiferromagnets, designed to capture the broad 
distribution of coupling constants in the effective model for the impurity degrees of freedom. Based on 
our results, we argue that in the presence of random magnetic couplings the conventional first-order 
spin-flop transition of an anisotropic antiferromagnet is split into two transitions at low temperatures, 
associated with separate order parameters along and perpendicular to the field axis. We demonstrate the 
existence of either a bicritcal point or a critical endpoint in the temperature-field phase diagram, with 
the consequence that signatures of the spin flop are more pronounced at elevated temperature. 



I. INTRODUCTION 

The interplay of magnetism and disorder is a fascinat- 
ing field of research in condensed matter physics. Various 
non-trivial low-temperature phases, like spin glasses, Bose 
glasses, random singlet phases etc., and associated phase 
transitions have been studied in both theory and experiment. 
A particularly interesting manifestation of quantum effects is 
impurity-induced magnetism in quantum paramagnets: The 
starting point is a Mott insulator with a finite spin gap which 
separates elementary spin excitations from the spin singlet 
ground state; examples are SrCu203, CuGeOs, PbNi2V208, 
SrCu2(B03)2, KCuCIs, or TICuCls. As demonstrated in 
recent experiments,-'"^ non-magnetic impurities, replacing 
magnetic ions in such quantum paramagnets, can induce 
effective magnetic moments. These impurity-induced mo- 
ments reveal themselves in a Curie-like behavior of the uni- 
form susceptibility, x ^ C/T, at intermediate temperatures. 
Remarkably, in the presence of three-dimensional couplings 
these induced moments can order at sufficiently low temper- 
atures, thus changing the spin-gapped paramagnetic ground 
state of the pure compound into an magnetically long-range 
ordered state upon doping. 

On the theoretical side, the appearance of effective mo- 
ments upon doping vacancies into the spin system is well 
understood in principle. It occurs in systems with confined 
spinons, i.e., elementary S = \ excitations; for spin 1/2 sys- 
tems it is best visualized in terms of broken singlet bonds 
where one spin is replaced by a vacancy. The liberated spin 
1/2 is confined to the vacancy at low energies, resulting in 
an effective spin 1/2 moment.^'-''' (In contrast, in host sys- 
tems with elementary S = 1/2 excitations, i.e., deconfined 
spinons, no moments are generated by introducing vacan- 
cies.) This theoretical picture has been supported by various 
numerical studies, in particular on spin chain-'^ and ladder 
systems. ^^"^^ 

A number of theoretical works also addressed impurity- 
induced antiferromagnetic ordering. On the one hand, nu- 
merical simulations on finite-size systems of spin gap mag- 
nets containing impurities studied signatures of magnetic 
ordering. ^^^'^ On the other hand, analytical approaches''®"^^ 
typically start from an effective spin-1/2 model involving the 



impurity-induced moments Si only: 

TLcs — '^^JijSi ■ Sj, (1) 

ij 

where the interaction Jij between two impurity moments 
depends on their distance as Jij oc exp(— rij/^), where 
^ is the magnetic correlation length of the host material. 
Due to the random locations of the impurities the system 
shows a broad distribution of coupling values Jij. Further- 
more, on bipartite lattices the sign of Jij will alternate as 
function of the Manhattan distance between i and j. This 
implies that classically all bonds can be satisfied with a Neel- 
type arrangement of the effective moments. In other words, 
Eq. (1) defines a strongly disordered non-frustrated quan- 
tum magnet. Real-space renormalization group studies^^'^-' 
indicate that the ground state of the model (1) with quan- 
tum spins 1/2 shows long-range order for any concentra- 
tion of impurity moments; this is supported by numerical 
simulations of vacancy-doped quantum paramagnets with 
confined spinons, which display magnetic order up to the 
percolation threshold. ■'^■'^ Thus, although the systems un- 
der consideration are long-range-ordered antiferromagnets, 
their properties can be expected to be strongly different 
from that of antiferromagnets without quenched disorder, 
due to the broadly distributed Jij. 

The purpose of this paper is to introduce a generalized 
mean-field theory, which takes into account the broad dis- 
tribution of coupling constants in Eq. (1). The central 
idea is to parameterize the spins according to their coupling 
strength to the environment, i.e., by their sum of coupling 
constants to all other spins, J, see Eq. (2) below. For each J 
a separate mean-field parameter will be introduced, leading 
to integral equations replacing the standard self-consistency 
relations. Although such a mean-field theory misses certain 
aspects of disorder physics, like localization phenomena, we 
will demonstrate that it captures various distinct properties 
of magnets described by Eq. (1). For instance, the overall 
behavior of the order parameter as function of tempera- 
ture is significantly different from non-disordered magnets 
and from conventional mean-field theory. In particular, we 
discuss the physics in an external field the presence of a 
magnetic anisotropy, relevant for most real materials. Here, 
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a spin-flop transition is expected to occur for fields par- 
allel to the easy axis, which indeed has been observed in 
TICui_a:Mga;Cl37 We argue that strong disorder leads to an 
interesting tennperature evolution of the spin-flop physics: 
at low temperature the transition is generically split into 
two (with at least one of them being continuous), whereas 
at elevated temperature a single first-order transition is re- 
stored. 

The bulk of the paper is organized as follows: In Sec. II 
we describe our generalized mean-field theory, together with 
the numerical procedure to solve the mean-field equations. 
Sec. Ill discusses the symmetries and possible phases of the 
model in the presence of a field parallel to the easy axis, and 
presents temperature-field phase diagrams obtained from 
the mean-field theory. In Sec. IV we take a closer look at 
the phase transitions, in particular at the spin-flop transi- 
tion. A comparison to experiments and available numerical 
results for vacancy-doped magnets is given in Sec. V. A 
brief outlook concludes the paper. 



II. MEAN-FIELD THEORY 

A. Parameterization 

In standard mean-field theory, the many-body problem 
is reduced to one single-spin problem in an effective field. 
Following this idea in the presence of disorder requires to 
consider distinct effective fields for all spins. Physicswise, 
we expect spins to behave differently if they have different 
couplings to their environment; spins in a similar environ- 
ment may behave in a similar fashion. This is the basis for 
the main simplification of our mean-field theory: We will 
parameterize the spins by the coupling sum J, defined as 



(2) 



Thus we replace the spin variables Si by S{J). The factor 
(— 1)*"-'"''^ accounts for the sublattice structure of the bi- 
partite lattice, with (i — j) being the Manhattan distance, 
and the second identity follows from the alternating sign 
of the coupling Jij. With this definition, Ji is the magni- 
tude of the effective field on spin i in a perfectly ordered 
antiferromagnetic state. 

For a description of the magnetic interactions in terms 
of the coupling sum, we introduce a probability distribution 
P(J) according to 



J^) 



(3) 



where N is the number of spins. Further, we need an inter- 
action function f{Ji,J2), defined as 
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NP{Ji)P{J2) 



Y,\jij\s{Ji-JiW2-jj 



(4) 



which is symmetric w.r.t. to Ji <->■ J2 and fulfills normal- 
ization conditions 

J dJ2P{J2)f{Jl,J2) = Jl. (5) 

With these definitions, a ferromagnetic Heisenberg 
model, H = — JijSi ■ Sj, takes the form 

W = -TVy" dJidJ2P{Ji)P{J2)f{Ji,J2)S{Ji) ■ S{J: 

Defining an effective field 

l(J)= j dJ2P{J2)f{J,J2){S{J2)), (7) 



(6) 



mi 



the mean-field Hamiltonian reads 



Hmf = j dJP{J) [to( J) + • S{J), (8) 

where we have included an external field B. 

For an antiferromagnet on a bipartite lattice different 
mean fields are required for the two sublattices A and B, 
which can become inequivalent in the presence of symme- 
try breaking and a finite field. Using Pa{J) = Pb{J) = 
P(J)/2, the effective field becomes 



mA,BiJ) = ± / dJ2P{J2)f{jj2) 



The corresponding mean-field Hamiltonian reads 



'Hml 



N 



dJP{J)( mA{J)+B Sa{J)+ 



mB(J) + B] Sb{J)) 



+ 



(9) 



(10) 



replacing (8). 

In general, two angles are needed to specify the orienta- 
tion of a spin. In the following, we will exclusively consider 
situations where H has a U(l) symmetry of rotations about 
the z axis (see Sec. III). Then we can describe the orien- 
tation of a spin just by one angle, (p, between spin and z 
axis; the angle in the xy plane can be fixed to zero, positive 
(negative) values of on the A (B) sublattice account for 
the staggered in-plane order. Hence, the expectation value 
of a spin can be written in the following way: 



{Sa,b{J))= {Sa,b{J)) 



sin 'PA,BiJ) 


cos (pA^BiJ) 



(11) 



The single-spin problem of Eq. (8) can be readily solved, 
yielding the mean-field equations: 



cos (Pa,b{J) = e"; 



[B + mA,B(J)] 
\B + mA,B{J)\ 



(12a) 
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and 



0a,b{J)) 



stanh 



B + mA,B{J) 



ksT/s 



(12b) 



The amplitude equation (12b) has been written for a spin 

with two states ±s, appropriate for quantum spins s = 
1/2 - this will be used in the numerical calculations. For 
continuous classical spins the tanh() needs to be replaced 
by a Brillouin function as usual. 

In the case of P{J) = 5{J — Jo) and /( Ji, J2) ~ Jo the 
mean-field equations (9,12) reduce to the self-consistency 
equation of the familiar Weiss mean-field theory, with 
Jo = zJ for a Hamiltonian with a nearest-neighbor cou- 
pling strength J and a coordination number of z. 



B. Magnetic anisotropy 

In this paper we consider the formally simplest source of 
magnetic anisotropy, namely an anisotropic exchange inter- 
action of easy-axis type (the behavior in the presence of a 
Dshyaloshinski-Moriya interaction is expected to be similar). 
Thus, in the Hamiltonian we perform the replacement 

JijSi . Sj - Jij [PiSfS^ + S^S]) + 5f 5|] (13) 

with an anisotropy constant /3 < 1; we keep the coupling 
to the external field as B ■ S. We note that in the context 
of impurity- induced magnetism, e.g., in TICuCIa some com- 
plications arise (which we will ignore here): The anisotropy 
of the effective interaction will depend on the external field 
and the interaction itself; furthermore the form of the field 
coupling will be modified, leading to an anisotropic g tensor. 

The anisotropy according to Eq. (13) requires the follow- 
ing replacement in the mean-field equations (12): 



mA,B 




(14) 



where m^^ and denote the components perpendicular 
and parallel to the easy axis. 



C. Choice of input parameters 

The described mean-field theory requires the coupling 
distribution P(J), Eq. (3), and the interaction function 
f{Ji,J2), Eq. (4), as input. Both functions are given by 
the underlying microscopic model. We have numerically 
determined P(J) and /(Ji, J2) for an effective model of 
impurity-induced order in quantum paramagnets, Eq. (1), 
as described in the Appendix. 

For the actual mean-field calculations we found it more 
convenient to use plausible model (i.e. fitting) functions in- 
stead. However, a difficulty arises here: P{J) and /(Ji, J2) 
cannot be chosen independently, because the normalization 
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FIG. 1: Model functions for the coupling-sum distribution 
P(J) obtained from different generating functions (;(Ji, J2): 

a) Lognormal-like distribution from Eq. (18), with the param- 
eters 0-2 = 0.699, fiL = -0.171, o-L = 1.3, and Jmax = 4.317. 

b) Gauss-like distribution from Eq. (19), with the parameters 
/til = 0.59, ai = 0.289, (T2 = 0.193, and Jmax = 1.929. In both 
cases, / dJP{J) = 1. 



conditions (5) cannot be easily fulfilled while preserving the 
symmetry w.r.t. J\ J2. We have therefore resorted 
to the following construction. From an arbitrary symmet- 
ric "generating" function g{Ji,J2) we define the functions 
P(J) and /(Ji, J2) according to: 



and 



f{Ji,J2) = a- 



(iJ2ff(J, J2) 



P(Jl)P(J2) 

Choosing the normalization factor as 



= / dJi 



rj 



dJ2giJi,J2), 



(15) 



(16) 



(17) 



all normalization conditions are fulfilled. 

As detailed in the Appendix, in vacancy-doped magnets 
the function P(J) will be peaked at a value which in- 
creases with impurity concentration, and it will be increas- 
ingly asymmetric at low concentrations, with a tail to higher 
J. Thus, among others, we have employed the following 
generating functions g{Ji, J2). 



^(Ji, J2) = exp ( - 



(In(Ji) - fiL? 



2al 



exp 



(ln(J2) - ml)- 
2ai 



X exp I — 



(J l - J2f 

2al 



d{Jmax — Jl)d{Jmax — J2) 

(18) 



generates a Lognormal-like distribution P(J) which cap- 
tures the coupling distribution for a small impurity concen- 
tration. In contrast, 



g{Ji,J2) = JiJ2exp 
(^1 



{■h-i>iY 



2al 



exp 



{■h - Ml) 
2<jI 



X exp 



2al 



0{Jmax — Jl)6{Jmax — J2) 

(19) 



leads to a Gauss-like distribution P(J), corresponding to 
higher impurity concentrations. The /i^ and fj; are free pa- 
rameters. Examples for the resulting P( J) are shown in Fig. 
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FIG. 2: (color online) Angle configurations 'fiiJ) occurring in 
the ordered phases of the antiferromagnet with field parallel to 
the easy axis, with blue/red (dark/gray) showing ^aI'-Pb- a) 
Ising phase, b) mixed phase, c) canted phase. The evolution 
fi-om a) to c) represents the behavior at low temperatures 
upon increasing the field, see Figs. 3,4 below. 



1; here and in the following the J values are scaled such that 
the mean value of ,/ is unity. In Fig. 8 below we show P( J) 
results of the numerical simulations for comparison. 

To model an easy-axis situation, we have chosen the 
anisotropy parameter /? = 0.9. 



III. PHASE DIAGRAMS 




A. Symmetries and phases 



FIG. 3: Temperature-field phase diagram of a disordered 
easy-axis antiferromagnet, obtained from the present mean- 
field theory, for the Lognormal-like distribution of coupling 
constants P(J) as shown in Fig. la. The calculation was 
done for spins with eigenvalues ±s, hence the axes are scaled 
by appropriate powers of s. First-order transitions (thick 
lines) were determined from crossing points of the free energy, 
whereas second-order transitions were obtained from extrap- 
olations of the order parameters. BP is a bi-critical point. 
Magnetization data along the dotted gray lines will be shown 
in Figs. 5, 7 below. 



We will restrict our attention to magnetic fields parallel 
to the easy axis in which case interesting spin-flop physics 
arises. Then, the symmetry of the antiferromagnetic Hamil- 
tonian (10) in the presence of a field is U(l) x Z2 (whereas 
for _B = we have U(l) x Z2 x Z2 which becomes SU(2) 
X Z2 in the absence of exchange anisotropy). Here we have 
assumed a unit cell size of 2 sites, and the Z2 symmetry 
corresponds to the exchange of the two sublattices. 

The possible phases are easily enumerated: (i) An Ising 
phase, where spins on the A (B) sublattice point preferren- 
tially up (down). This breaks the Z2 symmetry, but leaves 
the U(l) rotations about the z axis intact, (ii) A canted 
phase, where the spins order spontaneously perpendicular 
to the field, and have a finite component in the field direc- 
tion (equal for both sublattices) as well. Z2 and U(l) are 
broken, but a combination of sublattice exchange and 180° 
rotation is an intact Z2 symmetry, (iii) A mixed phase where 
Z2 and U(l) are fully broken, (iv) A disordered phase with 
no symmetry breaking. For non-zero field the spins point in 
the field direction only. 

The order parameters are the components of the stag- 
gered magnetization parallel and perpendicular to the direc- 
tion of the applied magnetic field, Afstagg,2 and A/stagg,xj/- 
Further, the ordered phases (i)-(iii) can be nicely character- 
ized by the angles Lp, Eq. (11). Without quenched disorder, 
i.e., for P(J) = (5(J — Jo), we have (i) Ising: Lp^ = 0, 
fB = 71"; (ii) canted: < (ps = —^A < tt; (iii) mixed: 
(fA < 0, (fB > 0, \(Pa\ 7^ IV-bI- Iri the disordered case 
with broad P{J), the angles become functions of J, with 
examples for pa,b{J) shown in Fig. 2. 



B. Numerical iteration of the mean-field equations 

For given functions P{J) and /(Ji, J2) and fixed values 
of applied magnetic field B, temperature T, and anisotropy 
constant /3 one can iterate the mean-field equations (9,12), 
using a linear discretization for the J values. The initial dis- 
tributions for the angles 'Pa,b{J) and amplitudes KS'a.b)! 
could be chosen random in principle; we found it more 
convenient to employ distributions corresponding to perfect 
Ising or XY order instead. The mean fields are calculated 
from Eq. (9); new amplitudes are obtained from Eq. (12b). 
Some care is required with the angle equation (12a) in the 
case of an Ising initial configuration: in addition to Eq. (12a) 
we used an update scheme where spins on the sublattice B 
flipping from ipB = tt to ips — are set to LpB = 7r/2 by 
hand in order to mix the symmetry breakings. 

These different initial conditions and different update 
schemes lead to potentially different fixed-point distribu- 
tions after convergence is reached. These correspond to 
different local minima in the free-energy landscape; com- 
parison of the free energies then yields the stable phase. 



C. Phase diagrams from mean-field theory 

In Figs. 3, 4 we show representative phase diagrams for 
the disordered easy-axis antiferromagnet with a longitudinal 
field, obtained from solving the mean-field equations (9,12). 
At zero field, the Ising order, present at low temperatures 
T, is destroyed at a continuous transition to a paramag- 
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FIG. 4: Temperature-field phase diagram as in Fig. 3, but for 
the Gauss-like distribution P(J) shown in Fig. lb. (CE is a 
critical endpoint.) 



netic phase. Applying a field to the Ising phase drives vari- 
ous transitions, resulting in a canted phase at intermediate 
fields and finally a field-polarized (disordered) phase at large 
field. The main difference to the text-book antiferromag- 
net is the presence of a mixed phase. The conventional 
first-order spin-flop transition is split into two transitions at 
low T: At some small field, there is a continuous transition 
from the Ising to a mixed phase, where an in-plane stag- 
gered magnetization Mstagg,xy perpendicular to the field 
develops. Only at a larger field, the Ising order measured by 
Mstagg,z is destroyed, leading to a spin-flop transition into 
a canted phase, see Figs. 5, 6. The origin of this behavior 
lies in the broadly distributed couplings: Already for small 
fields, weakly coupled spins with J <^ B cannot sustain the 
zero-field Ising order and flip in the field direction, while for 
strongly coupled spins (J ^ B) Ising order is favored. At 
low T, spins with J k B reach their lowest-energy state by 
canting - this results in an overall mixed phase and can be 
nicely seen in Fig. 2b. (However, at elevated temperatures 
the in-plane mean-field from these spins alone is not suffi- 
cient to establish a mixed phase.) At the spin-flop transition 
(mixed^canted at low T, Ising^canted at higher T), the 
spins with the largest couplings loose their Ising order. The 
evolution of the angle distributions fA.B{J) with increasing 
field at low T is shown in Fig. 2. 

Note that the part of the phase transition line between 
Ising and mixed phases at small T and B (dashed) is diffi- 
cult to extract numerically for a linear discretization of J. 
However, by a comparison of energies it is easy to prove 
that at T = 0, i? > the Ising phase is always unstable 
towards the mixed phase for distributions P(J) which are 
non-zero for arbitrarily small J. 

The two phase diagrams in Figs. 3, 4 differ in the exten- 
sions of the mixed phase, and in the character of the phase 
transition line between mixed and canted phase (see below). 
Clearly, the deviations from the textbook spin-flop behavior 
of non-disordered antiferromagnets are most pronounced for 
broad distributions P{J)- 
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FIG. 5: (color online) Magnetization curves for the 
Lognormal-like distribution P(J), Fig. la, at a) T/s'^ = 0.435 
and b) T/s^ = 2.177. The red (dashed) part indicates the 
jump at the first-order transition. 

Both the Neel temperature T/v and the critical field Be 
where long-range order is destroyed, scale with the mean 
value of the exchange constant (which is set to unity in 
our calculations). Similarly, the flop field Bfiop scales as 
(1-/3)1/2 times the mean exchange. Furthermore, an asym- 
metric distribution P{J) as the one in Fig. la leads to a 
larger T^, Be, Sfiop compared to a symmetric one with the 
same mean J, the reason being that the scales are primarily 
determined by the spins with large J. 



IV. PHASE TRANSITIONS AND CRITICAL 
PROPERTIES 

A. Magnetizations and critical exponents 

The Figs. 5, 6 show curves for the staggered magnetiza- 
tions (Afstagg.a;;^ and il/stagg.x) as well as the total magneti- 
zation M, along the gray dotted lines in the phase diagrams 
in Figs. 3, 4. The role of AIstagg,xy and Afstagg.z as order 
parameters can be clearly seen. 

Fig. 7 shows the total staggered magnetization vs. tem- 
perature for both coupling-sum distributions. The transi- 
tion to the disordered phase at high T (or B) is always 
continuous. The corresponding order-parameter exponent 
f3. defined through Afstagg oc (T^ - Tf, should be 1/2 
in conventional mean-field theory. Interestingly, Fig. 7a 
shows an overall behavior very different from this standard 
square-root law. This is again due to the broadly distributed 
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FIG. 6; (color online) As Fig. 5, but for the Gauss-like distri- 
bution P(J), Fig. lb, at a) T / = and b) T / = 0.675. 



couplings: Close to Tc the magnetization is effectively only 
carried by a small fraction of the spins with large J - note 
that the distribution P{J) corresponding to Fig. 7a has a 
pronounced tail at larger J. We note that asymptotically 
close to the phase transition standard mean-field behav- 
ior is restored within the present approach, with exponent 
/3=l/2. 



B. Spin-flop transition 

In a non-disordered easy-axis antiferromagnet, a first- 
order spin-flop transition from an Ising to a canted phase 
occurs upon increasing the field, leading to a jump in the 
total magnetization. 

In the present situation with quenched disorder, a mixed 
phase with both A/stagg,^^ and Afstagg.z non-zero occurs at 
low temperatures. The spin-flop transition splits; for small 
disorder (a narrow distribution P{J)) the mixed^canted 




FIG. 7: Staggered magnetization for a) the Lognormal-like 
distribution P{J), Fig. la, at B/s = 0.209, and b) the Gauss- 
like distribution P(J), Fig. lb, at B/s = 0.189. 



transition remains of first order, but becomes second-order 
at larger disorder, whereas the Ising-^mixed transition is al- 
ways continuous. As the mixed phase only exists at low T, 
the conventional first-order spin-flop transition is restored 
at elevated temperatures. This has the remarkable conse- 
quence that the jump in the magnetization is most pro- 
nounced at intermediate T, namely at the position of the 
critical end point or the bicritical point, respectively (see 
Figs. 5, 6). In general, the magnetization jump is larger in 
situations with less disorder because the collective spin flop 
is carried by a larger fraction of spins here. 



V. RELATION TO VACANCY-DOPED MAGNETS 

In the previous sections, we have described a general 
mean-field theory for disordered non-frustrated antiferro- 
magnets. We now discuss the applicability to impurity- 
induced magnetic order in quantum paramagnets. 



A. Experiments 

The low-energy physics of Mg-doped TICuCIs may be ex- 
pected to be well described by a Hamiltonian of the form 
(1), provided that the impurity concentration is small. Mag- 
netization measurements'^ have indicated the presence of a 
spin-flop transition at 1.8 K for fields parallel to the easy axis 
[2,0,1] at a field of approximately 0.35 T, significantly below 
the field corresponding to the bulk gap, 6T. Interestingly, 
the spin-flop field seems to be almost independent of the im- 
purity concentration (in the measured range of 0.8 - 2.5%). 
Our theory does not offer an easy explanation for that: As 
discussed above, the spin-flop field Sfiop is roughly propor- 
tional to the mean exchange (for fixed anisotropy parame- 
ter (3), which would result in a concentration dependence of 
Sfiop. Within the effective model (1), it is hard to envision 
a mechanism leading to a concentration-independent flop 
field, therefore we speculate that correlations between the 
impurities beyond this effective model play a role here. We 
also consider it possible that i?flop will actually decrease for 
impurity concentrations smaller than the ones measured. In 
this context we note that the experiments of Ref. 22 mapped 
out the phase diagram for TICuCIs doped with 1% Mg and 
showed the zero-field impurity-induced ordered phase to be 
continuously connected to the high-field bulk ordered phase. 
Together with recent theoretical studies^^'^^, this suggests 
that the impurity and bulk energy scales are not well sepa- 
rated there, i.e., the impurity concentration is too high to 
allow for a description using the effective model (1). 

Further magnetization measurements on Mg-doped 
TICuCIa would also be interesting regarding the temper- 
ature dependence of the spin-flop physics: Our theory pre- 
dicts interesting behavior for very low temperatures, where 
the spin-flop transition should split. This requires measure- 
ments down to e.g. 1/lOth of the ordering temperature T^; 
the experiments of Ref. 7 have T > T]\i/2. 
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B. Numerical results from Quantum Monte Carlo 
simulations 

Numerical simulations using Quantum Monte Carlo tech- 
niques can go beyond the effective model (1) and study 
the full system, i.e., quantum paramagnet plus vacancies. 
Those calculations have been reported in Refs. 14-17, but 
all were restricted to the spin-isotropic situation. These 
simulations mapped out the complete phase diagram, with 
distinct low-field and high-field ordered phases. Among the 
interesting aspects are the occurrence of Bose glass phases 
near the bulk field-ordered phase^^, and of a quantum dis- 
ordered phase at intermediate fields where impurity-induced 
transverse order is destroyed and the impurity moments ap- 
pear to form a random-singlet-like phase-^^. Clearly, these 
properties rely on the one hand on the quantum nature of 
the impurity-induced spin-1/2 moments and on the other 
hand on localization effects, both not captured by our mean- 
field approach. Quantum Monte Carlo calculations for the 
spin-anisotropic case studied by us would be interesting, but 
may be difficult due to the small energy scales involved in 
the spin-flop physics. 



VI. CONCLUSIONS 

We have proposed a mean-field theory for strongly disor- 
dered magnets, which takes into account the broad distri- 
bution of energy scales in the system. Parameterizing the 
spins by their sum of coupling constants, equivalent to the 
exchange field in a perfectly ordered state, yields a con- 
tinuous set of mean-field equations. We have applied the 
formalism to a model for impurity-induced order in spin-gap 
quantum magnets, and derived detailed phase diagrams as 
function of temperature and external field. We have shown 
that the conventional first-order spin-flop transition gener- 
ically splits at low temperatures, leaving room for a mixed 
phase with both transverse and longitudinal order. 

We envision our approach of continuous mean fields to 
be applicable to a number of interesting problems involving 
strong disorder, like magnetic ordering in dilute magnetic 
semiconductors^^, charge ordering in the presence of strong 
pinning, or electronic models treated within modifications 
of dynamical mean-field theory (DMFT)^^. 
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FIG. 8: Numerical results for the coupling-sum distribution 
P(J), in a host magnet with a correlation length of 7.9 lat- 
tice spacings, at impurity concentrations of a) 0.1% and b) 
1%. (The energy unit for the vertical axis is Jy of the bilayer 
model; note that the rather large values of J arise from the 
inaccurate short-range behavior of Eq. (A3).) 



APPENDIX A: COUPLING CONSTANTS FOR 
IMPURITY-INDUCED MOMENTS IN PARAMAGNETS 

Our mean-field theory requires the knowledge of the dis- 
tribution of the coupling-constant sums, P(J), and the in- 
teraction function f{Ji,J2)- We have numerically deter- 
mined these functions from an effective Hamiltonian for the 
impurity-induced moments of the form (1). 

The effective interaction between two impurity spins 
(sjj), coupled to two different spins (Sij) of a bulk system 
according to 

n = Hbuik + k(^S,- + S, ■ s,) (Al) 

can be determined in perturbation theory in K. In lowest 
order and in the static approximation, it is given by the 
to = bulk susceptibility: 

J^^ = K'{{Sr;Sf)){u, = 0). (A2) 

For vacancy-induced moments, K can be approximated by 

a bulk exchange coupling. 

To be specific, consider a bulk system consisting of dimers 
on a d-dimensional hypercubic lattice, with intra-dimer 
(inter-dimer) coupling J± (J||). Using a bond-operator for- 
malism in the linearized (harmonic) approximation^®, we 
find for the effective coupling in the long-distance limit: 

Jeff = - — {-l)C{-lY^(3-d)/22{d-l)/2^{d+l)/2 

2J|| (A3) 
X r-(<*-i)/2 exp(-r/^) 

where r denotes the distance between the impurity spins, 
(— 1)'' accounts the alternating sign of the interaction, and 

= (1) for spins i,j on the same (on different) sites 
of the dimer pairs. Note that the asymptotic behavior 
j--(d-i)/2 gxp(_^/^^ is generic, whereas the concrete value 
of ^ and the prefactor in (A3) depend on the level of ap- 
proximation used; in our case = J\\/{J± — 2(iJ||). For 
the purpose of our numerical simulation, we will employ 
Eq. (A3) with K = J|| for all distances r - this will over- 
estimate couplings at small r. The sign of Jgff leads to a 
non-frustrated system on a bipartite lattice. 

The simulation is done by randomly placing impurities on 
sites of a bilayer square lattice (i.e. d = 2), calculating P{J) 
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and /(Ji, J2), Eqs. (3,4), from the effective coupling con- 
stants (A3), and averaging the result over several impurity 
configurations. The bilayer system is assumed to be in its 
quantum disordered phase, close to the magnetic ordering 
transition. In the linearized bond-operator approach, the 
phase transition takes place at J\\/Ji. = 1/4; we choose in 
the following J\\/Ji_ = 0.249, corresponding to a correlation 



length of ^ « 7.9. Fig. 8 shows the resulting P{J) for two 
different impurity concentrations: In the low-concentration 
limit, the distribution is strongly asymmetric, with a tail to 
high energies. (Note that in this limit analytical results are 
available as well, see e.g. Ref. 23.) 
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